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Abstract. - Tlie "standard" Brownian motion master equation, used to describe thermal damping, 
is not completely positive, and does not admit a Monte Carlo method, important in numerical 
simulations. To eliminate both these problems one must add a term that generates additional 
position diffusion. He we show that one can obtain a completely positive simple quantum Brownian 
motion, efficiently solvable, without any extra diffusion. This is achieved by using a stochastic 
Schrodinger equation (SSE), closely analogous to Langevin's equation, that has no equivalent 
Markovian master equation. Considering a specific example, we show that this SSE is sensitive 
to nonlinearities in situations in which the master equation is not, and may therefore be a better 
model of damping for nonlinear systems. 
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All mechanical oscillators experience frictional damp- 
ing, in which they lose energy to their environment. This 
damping is accompanied by thermalization of the oscil- 
lator, since the environment is a large system, and thus 
a thermal bath. One would naturally like to model the 
effects of this damping on the oscillator, and indirectly 
on any other devices to which it is coupled, without hav- 
ing to describe the motion of the environment, with its 
many degrees of freedom. Frictional damping of quantum 
systems by an environment, often referred to as "quan- 
tum Brownian motion", has many applications. In par- 
ticular, it is essential in describing the behavior of nano- 
mechanical resonators, and is therefore of current inter- 
est in quantum nano-electro-mechanics (QNEMS) [1, 2] 
and quantum opto-mechanics [3,4]. For classical systems, 
Langevin's equation [5], containing only a deterministic 
frictional force and a Gaussian white noise source, gives 
a simple and excellent model of damping and thermal- 
ization. The situation is much less simple for quantum 
systems, however. 

The standard approach to obtaining a quantum equa- 
tion to model damping from a thermal bath is to derive a 
master equation for a linear oscillator by coupling it to a 
continuum of oscillators, and then trace out these oscilla- 
tors. This was first done by Caldeira and Leggett [6], and 
was brought to its completion in the tour-de-force by Hu 



et al. [7], building on the work of Unruh and Zurek [8]. 
They showed that an exact master equation for the oscil- 
lator alone could be derived, and that any non-markovian 
effects of the bath appear merely as time-dependent coef- 
ficients in this equation. This master equation, the exact 
BME, is [7] 

p = -iT{t) [x, {p, p}] ~ m [x, [x, p\] + Q{t) [x, b, p]]. (1) 

Here p is the density matrix for the oscillator, {p, p} = 
pp+pp is the anti-commutator, and x and p are the dimen- 
sionless position and momentum of the oscillator, defined 



X ^ {a + a^)/V2, 
p = -i{a-a^)/V2, 



(2) 
(3) 



where a is the oscillator's annihilation operator. The three 
functions of time, F, ^ and C, are the time-dependent coef- 
ficients that encode all non-Markovian effects of the bath, 
and are determined by the frequency of the resonator, w, 
the coupling to the bath (which sets the damping rate, 7), 
the structure of the bath (ohmic, super-ohmic, etc.), and 
the initial bath state (which is determined by the temper- 
ature, T). For reference, the physical position and mo- 
mentum are X = y^'h/{Lum)x and P = ^/Twornp, where m 
is the mass and lo the frequency of the oscillator. We note 
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that the master equation that results from using a bath 
of two-level systems, rather than oscillators, has the same 
form as Eq.(l) [9]. 

There are three issues with using the BME to model 
damping. The first is that the time-dependent coefficients 
are in general rather complex, and must be calculated nu- 
merically for a given problem [7] . The second is that this 
equation cannot be written in the Lindblad form [10]. Be- 
cause of this it has no equivalent stochastic Schrodinger 
equation (SSE, also referred to as an "unravelling") [11], 
and therefore cannot be simulated using the wave- function 
Monte Carlo method [12]. This method significantly re- 
duces the numerical overhead involved in simulating mas- 
ter equations, and is therefore of considerable practical 
importance. Thirdly, this master equation only generates 
the correct thermal steady-state for linear oscillators, and 
many questions of interest for mesoscopic oscillators in- 
volve non-lineairites of the oscillators themselves, or in- 
teractions with other nonlinear devices (such interactions 
effectively turn an oscillator into a nonlinear system). As 
far as the author is aware, no master equations have been 
derived for damping of nonlinear oscillators. Instead, the 
problem is avoided: since it is usually the case in practice 
that the nonlinearities are small compared to the harmonic 
motion, one assumes that the BME (and the SSE to be 
presented shortly) will give reasonable results. Ideally one 
would like to improve upon this situation, however, and 
with this in mind we consider the behavior of various BME 
variants, as well that of the SSE, for an example nonlinear 
system in detail below. 

The standard approach to addressing the first issue is to 
choose constant values for the time-dependent coefficients, 
so as to keep the resulting BME accurate to a good ap- 
proximation. Different choices for these values produce 
different variants. For oscillators at low temperatures, 
the relevant regime for QNEMS, a reasonable approximate 
version is [13] 



-i^[x,{P,P}]~ ^(2nT 



l)[x,[x,p]]. 



(4) 



Here ut is the average number of phonons in the resonator 
when it is at temperature T, being 



1 



nx 



eyip[huj/{kBT)] - 1 



(5) 



with /cb Boltzmann's constant. Note that in what follows, 
riT will always be defined as the specific function of T 
given in Eq.(5). We will refer to Eq.(4) as the "low tem- 
perature BME" (LBME). For completeness we note that a 
similar variant, sometimes refereed to in the quantum op- 
tics literature as the "standard" BME (SBME), is given 
by 

•7, 



-ij[xAp^p}] 



^ (!£)["' t^'^"- 



(6) 



This SBME is a reasonable approximation to the BME 
in the joint regime of high temperature and weak damp- 
ing [13]. These two variants still have no equivalent SSE. 



They also generate a non-physical term in the fluctuation 
spectrum of the resonator [14] (due to the fact that the 
evolution is not completely positive [15,16]), although this 
does not usually cause problems in the time domain. 

As first pointed out by Diosi [17], one can obtain a new 
BME that does have an equivalent SSE by adding a term 
to the LBME or SBME. The low-temperature version, 
which we will refer to here as the "completely positive" 
BME (PBME) is 



P 



-i-[x,{p,p}] 



16(2nT + l) 



{2nT + l)[x,[x,p]] 



_ 7 
4 

[pAp,p]] 



(7) 



The extra term generates diffusion in position. This equa- 
tion is a considerable improvement, as it is now a valid 
quantum evolution, and can be simulated using the wave- 
function Monte Carlo method. However, the PBME no- 
longer gives the correct thermal steady-state for a linear 
oscillator. This is because the steady-state average phonon 
number is increased by the position diffusion. Note that 
the position diffusion term in the PBME is not part of 
the exact BME derived by Hu et al. Nevertheless, a mas- 
ter equation of the form given in Eq.(7) can be derived 
as an approximation to the dynamics induced by a coUi- 
sional bath (rather than a linear bath of harmonic oscil- 
lators) [18-22,24]), and the high temperature version of 
Eq.(7) has also been derived as an approximate limit of a 
system interacting with a bath of oscillators [23] . 

We now present a quite different approach to mod- 
cling damping and thermal noise in mesoscopic oscilla- 
tors. This involves directly constructing a stochastic 
Schrodinger equation for this purpose, rather than first de- 
riving a master equation. In fact, while the resulting SSE 
is straightforward to simulate [25], it does not correspond 
to any Markovian master equation. The SSE satisfies all 
the properties required for the thermal damping of a lin- 
ear system, in particular generating the correct thermal 
steady-state. Since this equation has no additional posi- 
tion diffusion, it can be regarded as a genuine quantum 
version of simple damped Brownian motion. Further the 
motion of the mean values of position and momentum for 
a given trajectory are precisely those given by Langevin's 
equation (They contain a fiuctuating force that gener- 
ates momentum diffusion.) The SSE may therefore be 
regarded as a quantum equivalent of Langevin's equation. 
Performing simulations of a particular example, we show 
that while the SSE does not produce the correct (thermal) 
steady-state for nonlinear systems, it does better than the 
various versions of the Brownian motion master equation. 

The SSE we propose to model thermal damping is 



see Eq.(9). 



(8) 



Here H is the oscillator Hamiltonian, and the incre- 
ments dWi and dW2 are mutually independent Wiener 
noise increments satisfying the Ito calculus relation dW^ = 
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d|^) = 



H iriT I 2 2^ 



+ 2^7^ ( 2 + *^^^ ^ * 



IV') dt + 27 ( nT{x)x + nT{p)p + -7{p)x ) 1-0) di 



C 



dWi + 2y^jnT (- - iVpX + iC^pP 



dWo 



(9) 



dW2 = dt [11,26]. As is usual, the means, variances, and 
covariance of the position and momentum for the state jf/;) 
are denoted by {x) = {tp\x\^), {p) = {tplpltp), "14 = (a;^) — 
(x)2, Vp = (p2) - (p)2, and C^p = {xp + px)/2 - {x){p). 

To simplify the presentation, the above SSE is unnor- 
malized [11]. This means that to simulate the SSE, at 
each time-step dt one adds the increment dlip), and then 
must normalize the state before adding the next incre- 
ment. To obtain the density matrix for an oscillator un- 
dergoing Brownian motion, one simulates the above equa- 
tion for many realizations of the noise processes dWi . If 
we label the states resulting from each of the N noise pro- 
cess as |'0j(t)), j = 1, . . . , A^, then the density matrix for 
the system at time t is p{t) = (l/N) J2j lV'j(*))(V'j(*)l- 

The SSE (Eq. (9)) has the following properties: 

1. It gives a valid quantum evolution. 

2. It gives the Langevin equations for the mean values 
of position and momentum: 



d{x) 
dip) 



-'i{[x,H])dt, 
-i{[p,H])dt- 



-{p)dt + ^Jyn^dWi 



(10) 
(11) 



3. For a harmonic oscillator, the SSE reduces every 
initial state to a coherent state (with randomly fluctuating 
mean position and momentum) as i -^ 00. 

The first two properties confirm that the SSE realizes 
a quantum equivalent of classical damped Brownian mo- 
tion. The third shows that this evolution also gives the 
same final states as the standard Brownian motion master 
equation. Together these show us that the SSE provides 
a good description of simple damped quantum Brownian 
motion. The final property is now expected to follow, and 
we show below that it does: 

4. For a harmonic oscillator the SSE gives the expected 
thermal steady-state density matrix, being 



riT 



Vt)£ 



TlT 



riT + 1 



|n)(n| 



(12) 



where the states \n) are the phonon number states of the 
resonator. This can also be written as a Gaussian Wigner 
function. The steady-state covariances are 



X p cy 



TlT, 



(13) 



and Cll = 0. 

Constructing the SSE: To obtain an SSE with the four 
properties listed above, we begin by noting that for a linear 



oscillator, if we make two simultaneous continuous mea- 
surements, one of X and the other of p, this gives an SSE 
whose steady-state is a coherent state with randomly fluc- 
tuating mean position and momentum. The SSE, normal- 
ized this time, for this joint measurement is [11] 



see Eq.(15) 



(14) 



where k determines the rate at which the measure- 
ment extracts information (often called the "measurement 
strength" ) . The equations of motion for the means that 
result from this SSE are 



d{x) 
dip) 



?>kV^dWi 
ikVr,dW2 



SkC^pdW2, 



(16) 



%kCxpdWi, (17) 



= -^kV^ - ^kCip + 2fc, 


(18) 


= -SkV^ - SkClp + 2k, 


(19) 


= -8k{Va, + Vp)C,p. 


(20) 



and the equations for the covariances, when the wavefunc- 
tion for the state \ip) is Gaussian, are 



Vp 



For Brownian motion we require that the mean position 
is not driven by any noise, and that the mean momentum 
has the equation dip) = —{j/2)ip)dt + y^ynTdW2- We can 
produce the correct equations for the means by adding a 
linear feedback term, applied by a fictitious observer who 
has access to the stream of measurement results from the 
continuous measurements of x and p [27]. We choose the 
feedback Hamiltonian to cancel the noise driving the mean 
position, and to add damping to the mean momentum. 
The required feedback Hamiltonian is 



iJfb = 2hy/^^ Vp 



dW2 

~dr 



l{p) 



dWi 



dWi 
dt 



2h^^^[v^^ + a 



■^xp 



dW2 
dT 



P- 



With the addition of this fictitious feedback, there is only 
one more thing we need to do to obtain the Brownian 
motion SSE. This is to fix the value of the measurement 
strength, k. Since it is the thermal bath that causes 
both the measurement (decoherence) and the damping, 
k should be related to the damping rate, 7. To get this 
relationship right, we thus choose k so that the quadratic 
terms in our SSE (the terms that give the decoherence due 
to the measurement) match those of the simple forms of 
the standard Brownian motion master equation that give 
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d\i;) = -fc [{x - {x)f + {p-{p)f] dM) + V2fc [(.T - {x))dWi + {p-{p))dW2] IV) 



(15) 



the same damping rate [13,18]. (These quadratic terms 
also match the Lindblad master equation for the damp- 
ing of an optical cavity with the same damping rate [28]). 
This sets k = {'jnT)/2, and results in the SSE given by 
Eq.(9). 

By construction, the SSE gives a valid (completely posi- 
tive) quantum evolution (property 1 above) , and the equa- 
tions of motion for the expectation values of x and p satisfy 
property 2. It is well-established that continuous mea- 
surement of any linear combination of x and p produces 
Gaussian states as i -^ oo [29,30]. While to the authors 
knowledge no formal proof of this has been given to date, 
the proof is fairly straightforward, and follows from the 
results in [31]. Once the wave-function is Gaussian, the 
equations of motion for the covariances are those given in 
Eqs.(18)-(20), with the addition of the harmonic oscillator 
dynamics, given by 



K 



2luC^ 



= —2ujCxp, 
Cxp = L0{Vp~Vx). 



Vp 



(21) 
(22) 
(23) 



The the linear feedback Hamiltonian has no effect on the 
covariances. The resulting steady-state iaVx — Vp — 1/2, 
so the SSE satisfies property 4. 

To show the the SSE satisfies property 3, we first note 
that for a harmonic oscillator, in the long-time limit, the 
density matrix is a mixture of coherent states, due to prop- 
erty 4. The steady-state density matrix can therefore be 
obtained by calculating the probability density for (x) and 
(p) generated by Eqs.(lO) and (11). Solving these equa- 
tions is straightforward, and the resulting steady-state 
probability density for {x) and {p) is a Gaussian with zero 



mean and variances Vi 



{x) 



Vi 



ip) 



riT- The steady-state 



for the system is therefore a Gaussian mixture of coherent 
states, and has the variances Vx = Vp = 1/2 -I- ut- This is 
precisely the thermal state given by Eq.(12). 

As mentioned above, the SSE does not correspond to 
any Markovian master equation. The reason for this is 
the feedback. Without feedback one can average the SSE 
over many trajectories and obtain a Markovian master 
equation for the density matrix of the system. However, 
with feedback the dynamics of each trajectory depends on 
the state of the system during that trajectory in a way 
that cannot be reduced to a function of the density ma- 
trix (being the average over all trajectories) at the current 
time. The dynamics of the density matrix at time t cannot 
therefore be written purely in terms of the density matrix 
at t, making a description in terms of a Markovian master 
equation impossible. 

Damping of Nonlinear systems: We will now consider 
the behavior of the SSE for a particular example of a non- 
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Fig. 1: Here we show the evolution of the mean phonon num- 
ber of a harmonic oscillator, (n), under three evolutions that 
model damping and thermalization (Brownian motion). The 
temperature parameter in the equations is chosen so that the 
correct thermal steady-state value of (n) is 1. a) Solid line: the 
low temperature Brownian motion master equation (LBME); 
Dashed line: the completely positive version of the same mas- 
ter equation, b) Dark solid line: The stochastic Schrodinger 
equation described here, averaged over approximately 20,000 
trajectories; Light Solid Line: The result for the LBME, which 
is largely obscured behind. 



linear oscillator, and compare this with the behavior of the 
BME variants. As a reference we first solve the LBME and 
PBME master equations numerically, along with the SSE, 
for a linear oscillator. Measuring time in units of the oscil- 
lation period of the oscillator, we have uj = 2tt. We further 
choose the damping rate 7 = 4, the temperature so that 
riT = 1 (this means that T = A:B/(?iwln2)), and start the 
oscillator in the ground state (n = 0). In Fig. 1 we plot 
the time evolution of the mean phonon number under the 
two master equations and the SSE. The evolution of the 
SSE is indistinguishable from that of the LBME, both of 
which settle down to the expected steady-state value of 
the average phonon number, {n) = ut = 1. 

We now consider a x nonlinear oscillator, which has 
the Hamiltonian [33] 



_ff,ii = hLj{a'a) 



(24) 



This Hamiltonian has the same eigcnstates as the har- 
monic oscillator, \n), n = 0, 1,2,.. ., but the energy levels 
are now En = fiujn^ . It is only the two lowest energy levels 
that the x^^'' oscillator has in common with the harmonic 
oscillator (up to a shift of hijo/2). 

Simulating the nonlinear oscillator with the two master 
equations, with w, 7 and T (and thus ut) as before, one 
learns immediately that these give the same steady-state 
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Fig. 2: Here we show the steady-state mean energy, {E) oc (n^), 
for a X nonlinear oscillator, as well as the steady-state prob- 
ability distribution for the energy levels, n, for two evolutions; 
the low temperature Brownian motion master equation, and 
the SSE described here. The momentum damping rate is 7/2, 
and the temperature is chosen so that jit = 1, and a harmonic 
oscillator with the same lowest two energy levels has (n^) = 3. 
(a) circles: the value of (n^) given by the SSE; dark dashed 
line: thermal value for {n?); light dashed line: (n'^) = 3. (b)- 
(d): The distribution over the energy levels for three damping 
rates. Bars: the distribution given by the SSE; solid line: the 
thermal (Boltzmann) distribution; dashed line: the distribu- 
tion given by the Brownian motion master equations. 



distribution over the eigenstates, |n), as they do for the 
linear oscillator, and thus the same average value of n. 
That is, so long as the eigenstates are those of the har- 
monic oscillator, the steady-states of these master equa- 
tions know nothing about the energy levels of the system. 
This behavior is easily understood: the master equations 
eliminate the off-diagonal elements of the density matrix 
in the energy eigenbasis, and once this has happened the 
Hamiltonian has no effect on the evolution. Since the en- 
ergy eigenvalues only act on the dynamics via the Hamil- 
tonian, they can have no effect on the steady-state. 

The average energy for the nonlinear oscillator is no 
longer given by hu!{{n) + 1/2). The energy levels are 
now proportional to n^, and the average energy is (E) — 
huj{'n?). The steady-state distribution for the nonlinear 
oscillator given by the LBME is not the thermal (Boltz- 
mann) distribution. As noted above it remains the ther- 
mal distribution for the harmonic oscillator, and therefore 
gives too much weight to the higher energy levels. For 
T = k'Q/{'hLohi2) (that is, ut = 1), the thermal value of 
{n?) for the nonlinear oscillator is 0.49, and that given 
by the LBME is 3. To use the LBME to model damp- 
ing of the nonlinear oscillator at a given temperature, we 
could always alter the parameter ut so that the LBME 
gives the correct thermal value for (n^) at that temper- 
ature (and thus for the average energy). Of course this 
does not achieve the Boltzmann distribution. 

The behavior of the SSE for the nonlinear oscillator 
is quite different to that of the master equations. The 
SSE is sensitive to the energy levels of the oscillator, and 



the weighting of higher energy levels is suppressed as one 
would like. However, it also does not achieve the Boltz- 
mann distribution. In addition, the steady-state value of 
(n^), while usually closer to the thermal value than that 
given by the master equations, depends on the damping 
rate. In Fig. 2(a) we plot the steady-state value of {n?) 
as a function of the damping rate, with the parameter 
riT = 1 as before. For large damping the distribution given 
by the SSE tends to that given by the master equation 
(being the Boltzmann distribution for the harmonic oscil- 
lator). As the damping rate is reduced, {n?) decreases, 
and, after crossing the thermal value, settles into a limit- 
ing value, becoming independent of the damping rate for 
weak damping. For tit = 1, the SSE gives the thermal 
value for {-n?) when 7 ~ 0.8. For weak damping the SSE 
gives {v?) — 0.24 which is almost exactly half the thermal 
value. 

In Fig. 2 (b)-(d) we show the steady-state distribution 
over the energy levels given by the SSE for three values 
of the damping, and compare this to the correct thermal 
distribution for the nonlinear oscillator, and the thermal 
distribution for the harmonic oscillator. We see that for 
small damping, the SSE actually does the opposite of the 
LBME, in that it weights the lower energy levels more than 
the thermal distribution, and thus underestimates the av- 
erage energy. As the damping increases, the distribution 
changes and tends towards the thermal distribution for 
the harmonic oscillator, thus giving the same result as the 
LBME in the limit of large damping. We see from Fig. 2 
(c) that at no point does the distribution match the ther- 
mal distribution: when the average energy is equal to the 
thermal value, the SSE weights both n < 1 and n > 1 
more than the thermal distribution. 

The case of distant systems: It is worth noting that 
there is a wrinkle to using the SSE for simulating damp- 
ing on systems that are entangle with distant systems, 
and in which one wishes explicitly to examine questions of 
non-locality. In order to provide the feedback required to 
generate the damping, the bath (the fictitious observer) 
must have maximal knowledge of the state of the sys- 
tem. This is already implied by the fact that the SSE 
propagates a pure state, and causes no issue if the sys- 
tem is coupled only to other local systems. However, if 
one wishes to describe a thermal resonator that is entan- 
gled with a second, distant system, then the fact that we 
propagate a pure state for the resonator and the second 
system means that we give the bath full knowledge of this 
second system. This is also fine, unless we have an ob- 
server performing operations on the distant system while 
the resonator is evolving. In this case, allowing the bath 
access to the joint pure-state of both systems means that 
we are assuming infinitely fast communication from the 
distant observer to the local bath. 

If one does wish to use the SSE in such a situation, 
then one can preserve locality by modifying the simula- 
tion in the following way. One first includes the measure- 
ments of the distant observer in the SSE. Secondly, one 



K. Jacobs 



averages over these measurements to obtain a stochastic 
master equation (SME) that gives the state of knowledge 
of the bath. The feedback apphcd by the bath is now de- 
termined by the state-of-knowlcgde of the bath, and not 
by the SSE. The SSE and SME arc then propagated to- 
gether, both containing the feedback terms that arc now 
determined by the SME. This procedure preserves local- 
ity since the bath, and thus the feedback, only ever has 
access to local information. However, it should be noted 
that doing this does modify the feedback, and therefore 
the motion induced by the bath. That is, the motion 
is no-longer exactly the simple Brownian motion that we 
wish to simulate, at least while the statc-of-knowledgc of 
the bath remains appreciably mixed. The SSE is therefore 
not always appropriate for problems involving explicit in- 
vestigations of nonlocality. 

Conclusion: We have shown how to obtain a quantum 
model of simple damped Brownian motion that gives a 
completely positive evolution. This is achieved by con- 
structing a stochastic Schrodinger equation rather than a 
master equation, and by using the physical notion of feed- 
back from a measurement process, a concept borrowed 
from the subject of quantum feedback control [27]. 

While both the SSE and the Brownian motion master 
equations give the correct thermal distribution for the har- 
monic oscillator, interesting nano-mechnical devices are 
usually nonlinear, even if relatively weakly. We therefore 
explore using the LBME and the SSE to model damping 
of an oscillator in an extreme case, in which the Hamil- 
tonian is a x^^'' nonlinearity. While the LBME is com- 
pletely insensitive to the energy levels of this system (and 
thus, completely insensitive to the nonlinearity) the SSE 
is not. For weak damping the distribution given by the 
SSE is sensitive to the nonlinearity (and independent of 
the damping rate) although it overestimates the proba- 
bilities of the lower energy levels. In the limit of strong 
damping it reproduces the results of the LBME. For in- 
termediate damping the distribution of the energy levels, 
and thus the average energy, is dependent upon the damp- 
ing rate. These results suggest that the SSE is certainly 
no worse, and probably a better, model of thermal damp- 
ing for nonlinear systems than the various forms of the 
Brownian motion master equation. 

The present work highlights the fact that current tech- 
niques for modeling damping and thcrmalization in non- 
linear systems arc poor. It may be that stochastic 
Schrodinger equations, and the concept of feedback, will 
be a useful tool in developing more accurate models that 
are also efficient to implement. 
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